#ifndef VEC4D
#define VEC4D

#include <cmath>

struct vec4d {

  double x;
  double y;
  double z;
  double w;

  vec4d() : x(0.0), y(0.0), z(0.0), w(0.0) {
  }

  vec4d(double _x, double _y, double _z, double _w) : x(_x), y(_y), z(_z), w(_w) {
  }

  vec4d(const vec4d &b) : x(b.x), y(b.y), z(b.z), w(b.w) {
  }

  vec4d operator-(vec4d b) {
    return vec4d(x - b.x, y - b.y, z - b.z, w - b.w);
  }

  vec4d operator+(vec4d b) {
    return vec4d(x + b.x, y + b.y, z + b.z, w + b.w);
  }

  vec4d operator*(double value) const {
    return vec4d(x * value, y * value, z * value, w * value);
  }

  vec4d operator*(vec4d b) const {
    return vec4d(x * b.x, y * b.y, z * b.z, w * b.w);
  }

  vec4d operator/(double value) const {
    return vec4d(x / value, y / value, z / value, w / value);
  }

  vec4d operator/(vec4d b) const {
    return vec4d(x / b.x, y / b.y, z / b.z, w / b.w);
  }

  double &operator[](int index) {
    if (index == 0) return x;
    else if (index == 1) return y;
    else if (index == 2) return z;
    else return w;
  }

  bool operator==(vec4d b) {
    return (x == b.x && y == b.y && z == b.z && w == b.w);
  }

  bool operator!=(vec4d b) {
    return !(x == b.x && y == b.y && z == b.z && w == b.w);
  }

  double sqrLength(void) const { 
    return x * x + y * y + z * z + w * w; 
  }

  double length(void) const {
    return static_cast<double>(sqrt(this->sqrLength()));
  }

  double normalize(void) {
    double len =length();

    if (len > 0.0) {
      x = x / len;
      y = y / len;
      z = z / len;
      w = w / len;
    }
    return len;
  }
};

#endif